O objetivo desse notebook é aplicar a Linguagem R em Ciência de Dados na resolução de problemas de forma não supervisionada.
A base de dados USArrests contém estatísticas por 100.000 habitantes de prisões por agressão, assassinato e estupro em cada um dos 50 estados dos Estados Unidos em 1973.
data('USArrests')
usa = USArrests
head(usa)
É muito importante tratar valores ausêntes. A função is.na nos diz se um elemento é NA. Combinada com a função any, podemos descobrir se algum elemento da base é NA.
head(is.na(usa))
## Murder Assault UrbanPop Rape
## Alabama FALSE FALSE FALSE FALSE
## Alaska FALSE FALSE FALSE FALSE
## Arizona FALSE FALSE FALSE FALSE
## Arkansas FALSE FALSE FALSE FALSE
## California FALSE FALSE FALSE FALSE
## Colorado FALSE FALSE FALSE FALSE
print(any(is.na(usa)))
## [1] FALSE
Também é muito importante normalizar os dados. Mas poderiamos usar normalização por maxmin.
maxmin <- function(x) {
return (x - min(x))/(max(x) - min(x))
}
normalize <- function(df) {
return (sapply(df, maxmin))
}
usa = normalize(usa)
head(usa)
## Murder Assault UrbanPop Rape
## [1,] 12.4 191 26 13.9
## [2,] 9.2 218 16 37.2
## [3,] 7.3 249 48 23.7
## [4,] 8.0 145 18 12.2
## [5,] 8.2 231 59 33.3
## [6,] 7.1 159 46 31.4
Aqui vamos usar a normalização pela média e desvio padrão.
usa = scale(usa)
head(usa)
## Murder Assault UrbanPop Rape
## [1,] 1.24256408 0.7828393 -0.5209066 -0.003416473
## [2,] 0.50786248 1.1068225 -1.2117642 2.484202941
## [3,] 0.07163341 1.4788032 0.9989801 1.042878388
## [4,] 0.23234938 0.2308680 -1.0735927 -0.184916602
## [5,] 0.27826823 1.2628144 1.7589234 2.067820292
## [6,] 0.02571456 0.3988593 0.8608085 1.864967207
Agora nossos dados estão prontos. Vamos aplicar a função de agrupamento KMeans para 2 grupos.
cls = kmeans(usa, centers=2)
print(cls)
## K-means clustering with 2 clusters of sizes 20, 30
##
## Cluster means:
## Murder Assault UrbanPop Rape
## 1 1.004934 1.0138274 0.1975853 0.8469650
## 2 -0.669956 -0.6758849 -0.1317235 -0.5646433
##
## Clustering vector:
## [1] 1 1 1 2 1 1 2 2 1 1 2 2 1 2 2 2 2 1 2 1 2 1 2 1 1 2 2 1 2 2 1 1 1 2 2 2 2 2
## [39] 2 1 2 1 1 2 2 2 2 2 2 2
##
## Within cluster sum of squares by cluster:
## [1] 46.74796 56.11445
## (between_SS / total_SS = 47.5 %)
##
## Available components:
##
## [1] "cluster" "centers" "totss" "withinss" "tot.withinss"
## [6] "betweenss" "size" "iter" "ifault"
Analisando os grupos:
cls$cluster
## [1] 1 1 1 2 1 1 2 2 1 1 2 2 1 2 2 2 2 1 2 1 2 1 2 1 1 2 2 1 2 2 1 1 1 2 2 2 2 2
## [39] 2 1 2 1 1 2 2 2 2 2 2 2
Conseguimso descobrir qual o número de grupos ideal? Uma alternativa é utilizar uma medida de validação de agrupamento como a silhueta. Para isso precisamos carregar a biblioteca cluster.
require('cluster')
## Loading required package: cluster
sil = silhouette(cls$cluster, dist(usa))
print(sil)
## cluster neighbor sil_width
## [1,] 1 2 0.3354254
## [2,] 1 2 0.3254396
## [3,] 1 2 0.3911248
## [4,] 2 1 0.1345762
## [5,] 1 2 0.3774168
## [6,] 1 2 0.2972633
## [7,] 2 1 0.5306315
## [8,] 2 1 0.1743774
## [9,] 1 2 0.5035739
## [10,] 1 2 0.4044931
## [11,] 2 1 0.4194608
## [12,] 2 1 0.5680837
## [13,] 1 2 0.3354701
## [14,] 2 1 0.4428995
## [15,] 2 1 0.5945178
## [16,] 2 1 0.5362530
## [17,] 2 1 0.3496924
## [18,] 1 2 0.4387476
## [19,] 2 1 0.5625723
## [20,] 1 2 0.4890198
## [21,] 2 1 0.3837284
## [22,] 1 2 0.5165960
## [23,] 2 1 0.5986325
## [24,] 1 2 0.2620747
## [25,] 1 2 0.1134541
## [26,] 2 1 0.5240161
## [27,] 2 1 0.5926288
## [28,] 1 2 0.4299099
## [29,] 2 1 0.5870108
## [30,] 2 1 0.1782450
## [31,] 1 2 0.5222275
## [32,] 1 2 0.3865530
## [33,] 1 2 0.2596933
## [34,] 2 1 0.5143664
## [35,] 2 1 0.3618181
## [36,] 2 1 0.4064074
## [37,] 2 1 0.1942476
## [38,] 2 1 0.5253841
## [39,] 2 1 0.3635329
## [40,] 1 2 0.3650274
## [41,] 2 1 0.5325247
## [42,] 1 2 0.3175580
## [43,] 1 2 0.3471312
## [44,] 2 1 0.4078063
## [45,] 2 1 0.4408083
## [46,] 2 1 0.2622062
## [47,] 2 1 0.3313743
## [48,] 2 1 0.4592799
## [49,] 2 1 0.5891363
## [50,] 2 1 0.4400334
## attr(,"Ordered")
## [1] FALSE
## attr(,"call")
## silhouette.default(x = cls$cluster, dist = dist(usa))
## attr(,"class")
## [1] "silhouette"
Podemos mostrar graficamente o valor da silhueta por amostra.
plot(sil)
Ou usar a média para avaliar se o valor é bom.
mean(sil[,3])
## [1] 0.408489
Podemos variar o número de grupos para tentar encontrar o valor ótimo de grupos baseado na medida silhueta.
cls = sapply(2:10, function(k) {
aux = kmeans(usa, centers=k)
sil = silhouette(aux$cluster, dist(usa))
mean(sil[,3])
})
Podemos plotar os valores da silhueta.
names(cls) = 2:10
plot(names(cls), cls, type = "l")
Usando a biblioteca plotly.
require(plotly)
## Loading required package: plotly
## Loading required package: ggplot2
##
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
##
## last_plot
## The following object is masked from 'package:stats':
##
## filter
## The following object is masked from 'package:graphics':
##
## layout
data = data.frame(cls, k=2:10)
plot_ly(data, x=~k, y=~cls, type='scatter', mode='lines')
Vamos carregar duas bibliotecas para nos ajudar a visualizar os dados.
require('gridExtra')
## Loading required package: gridExtra
require('factoextra')
## Loading required package: factoextra
## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBa
Vamos visualizar os agrupamentos.
plot = lapply(2:5, function(k) {
cls = kmeans(usa, centers=k)
fviz_cluster(cls, data=usa)
})
do.call("grid.arrange", plot)
Um experimento interessante é a tarefa de comprimir áudio utilizando o algoritmo de agrupamento KMeans. Para isso vamos precisar das bibliotecas audio e seewave.
require('audio')
## Loading required package: audio
require('seewave')
## Loading required package: seewave
##
## Attaching package: 'seewave'
## The following object is masked from 'package:plotly':
##
## export
Vamos baixar uma música da internet com permissão de reprodução pelo pixabay. Essa música tem o formato MP3.
url = 'https://cdn.pixabay.com/download/audio/2021/04/07/audio_0bed53371b.mp3?filename=blues-vibes-100-bpm-michael-kobrin-3780.mp3'
download.file(url, 'audio.mp3')
Na sequencia, precisamos converter para WAV e pegar uma amostra para podermos demostrar a compressão. Para isso vamos usar o programa ffmpeg.
system('ffmpeg -i audio.mp3 audio.wav')
system('ffmpeg -i audio.wav -ac 1 mono.wav')
system('ffmpeg -i mono.wav -ss 00:00 -to 00:10 out.wav')
Podemos carregar a música.
sound = load.wave("out.wav")
class(sound)
## [1] "audioSample"
E plotar…
plot(sound, type = "l")
Usando plotly…
data = data.frame(as.numeric(sound), x=1:length(sound))
plot_ly(data, x=~x, y=~sound, type='scatter', mode='lines')